1. Field of the Invention
The present invention relates to an image processing apparatus and a computer program product that estimate an amount of dye in a stained sample from a stained sample image that is obtained by imaging the stained sample stained with a plurality of dyes.
2. Description of the Related Art
As one of physical quantities that indicate inherent physical property of a subject of imaging, there is a spectral transmittance. The spectral transmittance is a physical quantity that indicates a rate of transmitted light of incident light at each wavelength, and is inherent information for a substance, and the value thereof is not affected by an extrinsic influence, unlike color information that is dependent on a change of illumination light, as an RGB value. Therefore, the spectral transmittance is used in various fields as information to reproduce color of the subject. For example, for a living tissue specimen, particularly in a field of pathology using pathological samples, a technology of estimating spectral transmittance is used in the analysis of images of samples as one example of spectroscopic characteristics.
In pathology, such a process is widely practiced that a pathological sample is magnified to be observed using a microscope after slicing a block sample obtained by excision of an organ or a sample obtained by needle biopsy into piece having several microns of thickness to obtain various findings. Transmission observation using an optical microscope is one of the methods that are most widely used because the equipments are relatively inexpensive and easy to be handled, and this method has a long history. In this case, because a sliced sample absorbs or scatters little light and is almost transparent and colorless, it is common to stain the sample with a dye prior to observation.
As staining methods, various methods have been proposed, and there are more than 100 methods in total. Particularly for pathological samples, hematoxylineosin stain (hereinafter, “H&E stain”) using bluish purple hematoxylin and red eosin as dyes is generally used.
Hematoxylin is a natural substance that is extracted from a plant, and has no stainability itself. However, hematin, which is an oxide of hematoxylin, is a basophilic dye and combines with a substance negatively charged. Because deoxyribonucleic acid (DNA) included in a cell nucleus is negatively charged due to a phosphate group included therein as a structural element, DNA combines with hematin to be stained bluish purple. As described, substance having stainability is not hematoxylin but its oxide, hematin; however, because it is common to use hematoxylin as the name of the dye, this practice is followed in the following explanation. On the other hand, eosin is acidophilic dye, and combines with a substance positively charged. Amino acid and protein are charged positively or negatively depending on a pH environment, and have inclination to be charged positively under acidity. For this reason, there is a case where acetic acid is added to eosin. Protein included in a cytoplasm combines with eosin to be stained red or light red.
In a sample subjected to H&E stain (stained sample), a cell nucleus, bone tissues, and the like are stained bluish purple, and cytoplasm, a connective tissue, red corpuscles, and the like are stained red, to become easily visible. As a result, an observer can grasp the size, positional relation, or the like of elements structuring a cell nuclei or the like, and can determine a state of the sample morphologically.
Observation of samples is performed by multiband imaging the sample to be displayed on a display screen of an external device, other than visual inspection by an observer. In the case of displaying on a display screen, processing to estimate spectral transmittance at each sample point from the multiband image obtained, processing to estimate an amount of dye with which the sample is stained based on the estimated spectral transmittance, processing to correct color of the image based on the estimated amount of dye, and the like are performed. As a result, variation in a property of the camera, a stain condition, or the like are corrected, and an RGB image for display of the sample is generated. FIG. 19 is a view showing one example of a composed RGB image. If the estimation of an amount of dye is appropriately performed, a sample that is stained dark or the sample stained light can be corrected to an image in colors equivalent to the sample that is properly stained.
As a method of estimating spectral transmittance at each sample point from multiband images of the samples, for example, an estimation method by principal component analysis (for example, “Development of support systems for pathology using spectral transmittance—The quantification method of stain conditions”, Proceedings of SPIE, Vol. 4684, 2002, pp. 1516-1523), an estimation method by Wiener's estimation (for example, “Color Correction of Pathological Images Based on Dye Amount Quantification”, OPTICAL REVIEW, Vol. 12, No. 4, 2005, pp. 293-300), and the like can be used. Wiener's estimation is widely known as one of liner filtering methods by which an original signal is estimated from an observed signal on which noise is superimposed, and is a method in which minimization of error is performed considering statistical properties of a subject of observation and characteristics of noise (observation noise). Because some noise is included in a signal from a camera, Wiener's estimation is a very effective as a method of estimating an original signal.
A method of synthesizing an RGB image from a multiband image of a sample is explained. First, a multiband image of a sample is obtained. For example, using a technique disclosed in Japanese Patent Laid-Open Publication No. H7-120324, multiband images are taken by a frame sequential method while switching 16 pieces of band-pass filters by rotating a filter wheel. Thus, multiband images having pixel values of 16 bands at each sample point can be obtained. Although a dye is three-dimensionally distributed in the sample being a subject of observation in an actual state, it cannot be taken as a three-dimensional image as it is with an ordinary transmission observing system, and is observed as a two-dimensional image in which illumination light that has passed the sample is projected on an imaging device of the camera. Accordingly, each point mentioned herein signifies a point on the sample corresponding to each projected pixel of the imaging device.
For an arbitrary point x in the imaged multiband image, there is relation expressed as in the following equation (1) based on a response system of the camera, between a pixel value g(x, b) in band b and spectral transmittance t(x, λ) of a corresponding point on the sample.g(x,b)=∫λf(b,λ)s(λ)e(λ)t(x,λ)dλ+n(b)  (1)where λ indicates wavelength, f(b, λ) indicates spectral transmittance of a b-th filter, s(λ) indicates spectral sensitivity property of the camera, e(λ) indicates spectral radiance property of illumination light, and n(b) indicates observation noise in band b. The variable b is a serial number to identify a band, and is an integer that satisfies 1≦b≦16 in this example.
In an actual calculation, the following equation (2) obtained by discretizing equation (1) is used.G(x)=FSET(x)+N  (2)
When the number of sample points in a direction of wavelength is D and the number of bands is B (B=16 in this example), G(x) is a matrix of B×1 corresponding to the pixel value g(x, b) at the point x. Similarly, T(x) is a matrix of D×1 corresponding to t(x, λ), and F is a matrix of B×D corresponding to f(b, λ). On the other hand, S is a diagonal matrix of D×D and a diagonal element corresponds to s(λ). Similarly, E is a diagonal matrix of D×D and a diagonal element corresponds to c(λ). N is a matrix of B×1 corresponding to n(b). In equation (2), because expressions of a plurality of bands are put together using a matrix, a variable b indicating a band is not specified in equation (2). Moreover, an integral of the wavelength λ is replaced with the product of matrices.
To simplify description, a matrix H defined by the following equation (3) is introduced. H is also called a system matrix.H=FSE  (3)
Next, spectral transmittance at each sample point is estimated from the imaged multiband image using Wiener's estimation. An estimation value of spectral transmittance (spectral transmittance data) {circumflex over (T)}(x) can be calculated by the following equation (4).{circumflex over (T)}(x)=WG(x)  (4)where W is expressed by a following equation, and is called “Wiener's estimation matrix” or “estimation operator used in Wiener's estimation”. In the explanation below, W is simply referred to as “estimation operator”.W=RSSHt(HRSSHt+RNN)−1  (5)where ( )t indicates a transposed matrix, and ( )−t indicates an inverse matrix. Furthermore, RSS is a matrix of D×D, and expresses an autocorrelation matrix of spectral transmittance of the sample. RNN is a matrix of B×B, and expresses an autocorrelation matrix of noise of the camera used for imaging.
After thus estimating the spectral transmittance data {circumflex over (T)}(x), an amount of dye at a corresponding sample point (corresponding point) is estimated based on this {circumflex over (T)}(x). Dyes to be subjects of estimation are three kinds of dyes: hematoxylin, eosin that stains cytoplasm, and eosin that stains red corpuscles or red corpuscles that are not stained. Three kinds of dyes are abbreviated as dye H, dye E, and dye R, respectively. Precisely, red corpuscles have a peculiar color even in a not stained state, and after the H&E stain is performed, the color of red corpuscles and the color of eosin that has changed in a staining process are superimposed with each other at the time of observation. Therefore, in precise, color obtained by combining the both is called dye R.
Generally, it is known that Lambert-Beer law expressed by the following equation (6) is satisfied between intensity I0(λ) of incident light at each wave length λ and intensity I(λ) of emitting light in a substance passing light.
                                          I            ⁢                                                  ⁢                          (              λ              )                                                          I              0                        ⁡                          (              λ              )                                      =                  ⅇ                                    -                              k                ⁡                                  (                  λ                  )                                                      ·            d                                              (        6        )            where k(λ) indicates an inherent value of a substance dependent on wavelength, and d indicates thickness of a substance. Moreover, the left side of equation (6) indicates spectral transmittance.
When a sample subjected to H&E stain is stained with three kinds of dyes of dye H, dye E, and dye R, a following equation (7) is satisfied at each wavelength λ by Lambert-Beer law.
                                          I            ⁡                          (              λ              )                                                          I              0                        ⁡                          (              λ              )                                      =                  ⅇ                      -                          (                                                                                          k                      H                                        ⁡                                          (                      λ                      )                                                        ·                                      d                    H                                                  +                                                                            k                      E                                        ⁡                                          (                      λ                      )                                                        ·                                      d                    E                                                  +                                                                            k                      R                                        ⁡                                          (                      λ                      )                                                        ·                                      d                    R                                                              )                                                          (        7        )            where kH(λ), kE(λ), and kR(λ) indicate k(λ) corresponding to dye H, dye E, and dye R, respectively, and are standard spectral properties of respective dyes that stain the sample, for example. Furthermore, dH, dE, and dR indicate virtual thickness of dye H, dye E, and dye R at each sample point corresponding to each image position of the multiband image. Originally, dyes are dispersed in a sample, and therefore, thickness is not a correct idea. However, this can be an index of a relative amount of dye that indicates how much amount of dye is present compared to a case where the sample is stained with a single dye. In other words, dH, dE, and dR indicate amounts of dye H, dye E, and dye R, respectively. The values kH(λ), kE(λ), and kR(λ) can be easily acquired from Lambert-Beer law, by preparing samples that are stained respectively using dye H, dye E, and dye R, and by measuring spectral transmittance with a spectrometer.
If logarithms of both sides of equation (7) are taken, the flowing equation (8) is obtained.
                                          -            log                    ⁢                                          ⁢                                    I              ⁡                              (                λ                )                                                                    I                0                            ⁡                              (                λ                )                                                    =                                                            k                H                            ⁡                              (                λ                )                                      ·                          d              H                                +                                                    k                E                            ⁡                              (                λ                )                                      ·                          d              E                                +                                                    k                R                            ⁡                              (                λ                )                                      ·                          d              R                                                          (        8        )            
When an element corresponding to the wavelength λ of the spectral transmittance data {circumflex over (T)}(x) thus estimated is {circumflex over (t)}(x, λ), and if this is substituted in equation (8), the flowing equation (9) is obtained.−log {circumflex over (t)}(x,λ)=kH(λ)·dH+kE(λ)·dE+kR(λ)·dR  (9)
Estimated absorbance â(x, λ) can be calculated according to the following equation (10) based on the spectral transmittance {circumflex over (t)}(x, λ).â(x,λ)=−log {circumflex over (t)}(x,λ)  (10)
Therefore, equation (9) can be replaced with the following equation (11).â(x,λ)=kH(λ)·dH+kE(λ)·dE+kR(λ)·dR  (11)
In equation (11), unknown variables are three variables of dH, dE and dR. Therefore, if simultaneous equations are acquired from equation (11) for at least three different wavelengths λ, these can be solved. To further improve the accuracy, acquiring simultaneous equations from equation (11) for four or more different wavelengths λ, multiple regression analysis can be performed. For example, simultaneous equations acquired from equation (11) for three wavelengths λ1, λ2, λ3 can be expressed in a matrix as the following equation (12).
                              (                                                                                          a                    ^                                    ⁡                                      (                                          x                      ,                                              λ                        1                                                              )                                                                                                                                            a                    ^                                    ⁡                                      (                                          x                      ,                                              λ                        2                                                              )                                                                                                                                            a                    ^                                    ⁡                                      (                                          x                      ,                                              λ                        3                                                              )                                                                                )                =                              (                                                                                                      k                      H                                        ⁡                                          (                                              λ                        1                                            )                                                                                                                                  k                      E                                        ⁡                                          (                                              λ                        1                                            )                                                                                                                                  k                      R                                        ⁡                                          (                                              λ                        1                                            )                                                                                                                                                              k                      H                                        ⁡                                          (                                              λ                        2                                            )                                                                                                                                  k                      E                                        ⁡                                          (                                              λ                        2                                            )                                                                                                                                  k                      R                                        ⁡                                          (                                              λ                        2                                            )                                                                                                                                                              k                      H                                        ⁡                                          (                                              λ                        3                                            )                                                                                                                                  k                      E                                        ⁡                                          (                                              λ                        3                                            )                                                                                                                                  k                      R                                        ⁡                                          (                                              λ                        3                                            )                                                                                            )                    ⁢                      (                                                                                d                    H                                                                                                                    d                    E                                                                                                                    d                    R                                                                        )                                              (        12        )            
Equation (12) is replaced with the following equation (13).{circumflex over (A)}(x)=Kd(x)+ε  (13)
When the number of sample points in a direction of wavelength is D, Â(x) is a matrix of D×1 corresponding to â(x, λ), K is a matrix of D×3 corresponding to k(λ), d(x) is a matrix of 3×1 corresponding to dH, dE, and dR at the point x, and ε is a matrix of D×1 corresponding to an error.
According to equation (13), the amount of dye dH, dE, and dR are calculated using a least square method. The least square method is a method of determining d(x) such that the square sum of the error is minimized in single regression analysis, and it can be calculated by the following equation (14).d(x)=(KTK)−1KTÂ(x)  (14)
If the amount of dye dH, dE, and dR are acquired as described above, a change in the amount of dye in the sample can be simulated by correcting these amounts. Specifically, it is adjusted by multiplying the respective amount of dye dH, dE, and dR by appropriate coefficients αH, αE, αR, to be substituted in equation (7). Thus, new spectral transmittance t*(x, y) can be obtained by the following equation (15).t*(x,λ)=e−(kH(λ)·αHdH+kE(λ)·αEdE+kR(λ)·αRdR)  (15)
If equation (15) is substituted in the equation (1), an image of the sample in which amount of dye is virtually changed can be synthesized. In this case, it can be calculated assuming noise n(b) is zero.
By estimating the amount of dye at the arbitrary point x in a multiband image by the above procedure, the amount of dye of the sample can be corrected by virtually adjusting the amount of dye at each sample point and synthesizing an image of the sample after adjustment. Therefore, even if there is variation in stain of the sample, for example, a user can observe an image that is adjusted to an appropriate stain condition.